Load all required libraries.

library(tidyverse)
## Warning: package 'tidyverse' was built under R version 3.6.3
## -- Attaching packages -------------------------------------------------------------------------- tidyverse 1.3.0 --
## v ggplot2 3.3.2     v purrr   0.3.4
## v tibble  3.0.3     v dplyr   1.0.0
## v tidyr   1.1.0     v stringr 1.4.0
## v readr   1.3.1     v forcats 0.5.0
## Warning: package 'ggplot2' was built under R version 3.6.3
## Warning: package 'tibble' was built under R version 3.6.3
## Warning: package 'readr' was built under R version 3.6.3
## Warning: package 'dplyr' was built under R version 3.6.3
## Warning: package 'forcats' was built under R version 3.6.3
## -- Conflicts ----------------------------------------------------------------------------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(plotly)
## Warning: package 'plotly' was built under R version 3.6.3
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(broom)
## Warning: package 'broom' was built under R version 3.6.3

Read in raw data from RDS.

raw_data <- readRDS("./n1_n2_cleaned_cases.rds")

Make a few small modifications to names and data for visualizations.

final_data <- raw_data %>% mutate(log_copy_per_L = log10(mean_copy_num_L)) %>%
  rename(Facility = wrf) %>%
  mutate(Facility = recode(Facility, 
                           "NO" = "WRF A",
                           "MI" = "WRF B",
                           "CC" = "WRF C"))

Seperate the data by gene target to ease layering in the final plot

#make three data layers
only_positives <<- subset(final_data, (!is.na(final_data$Facility)))
only_n1 <- subset(only_positives, target == "N1")
only_n2 <- subset(only_positives, target == "N2")
only_background <<-final_data %>% 
  select(c(date, cases_cum_clarke, new_cases_clarke, X7_day_ave_clarke, cases_per_100000_clarke)) %>%
  group_by(date) %>% summarise_if(is.numeric, mean)

#specify fun colors
background_color <- "#7570B3"
seven_day_ave_color <- "#E6AB02"
marker_colors <- c("N1" = '#1B9E77',"N2" ='#D95F02')
#remove facilty C for now
#only_n1 <- only_n1[!(only_n1$Facility == "WRF C"),]
#only_n2 <- only_n2[!(only_n2$Facility == "WRF C"),]

only_n1 <- only_n1[!(only_n1$Facility == "WRF A" & only_n1$date == "2020-11-02"), ]
only_n2 <- only_n2[!(only_n2$Facility == "WRF A" & only_n2$date == "2020-11-02"), ]

Build the main plot

      #first layer is the background epidemic curve
        p1 <- only_background %>%
              plotly::plot_ly() %>%
              plotly::add_trace(x = ~date, y = ~new_cases_clarke, 
                                type = "bar", 
                                hoverinfo = "text",
                                text = ~paste('</br> Date: ', date,
                                                     '</br> Daily Cases: ', new_cases_clarke),
                                alpha = 0.5,
                                name = "Daily Reported Cases",
                                color = background_color,
                                colors = background_color,
                                showlegend = FALSE) %>%
            layout(yaxis = list(title = "Clarke County Daily Cases", showline=TRUE)) %>%
            layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
        
        #renders the main plot layer two as seven day moving average
        p1 <- p1 %>% plotly::add_trace(x = ~date, y = ~X7_day_ave_clarke, 
                             type = "scatter",
                             mode = "lines",
                             hoverinfo = "text",
                            text = ~paste('</br> Date: ', date,
                                                     '</br> Seven-Day Moving Average: ', X7_day_ave_clarke),
                             name = "Seven Day Moving Average Athens",
                             line = list(color = seven_day_ave_color),
                             showlegend = FALSE)
      

        
        #renders the main plot layer three as positive target hits
        
        p2 <- plotly::plot_ly() %>%
          plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
                                       type = "scatter",
                                       mode = "markers",
                                       hoverinfo = "text",
                                       text = ~paste('</br> Date: ', date,
                                                     '</br> Facility: ', Facility,
                                                     '</br> Target: ', target,
                                                     '</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
                                       data = only_n1,
                                       symbol = ~Facility,
                                       marker = list(color = '#1B9E77', size = 8, opacity = 0.65),
                                       showlegend = FALSE) %>%
          plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
                                       type = "scatter",
                                       mode = "markers",
                                       hoverinfo = "text",
                                       text = ~paste('</br> Date: ', date,
                                                     '</br> Facility: ', Facility,
                                                     '</br> Target: ', target,
                                                     '</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
                                       data = only_n2,
                                       symbol = ~Facility,
                                       marker = list(color = '#D95F02', size = 8, opacity = 0.65),
                                       showlegend = FALSE) %>%
            layout(yaxis = list(title = "SARS CoV-2 Copies/L", 
                                 showline = TRUE,
                                 type = "log",
                                 dtick = 1,
                                 automargin = TRUE)) %>%
            layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
        
        #adds the limit of detection dashed line
        p2 <- p2 %>% plotly::add_segments(x = as.Date("2020-03-14"), 
                                          xend = ~max(date + 10), 
                                          y = 3571.429, yend = 3571.429,
                                          opacity = 0.35,
                                          line = list(color = "black", dash = "dash")) %>%
          layout(annotations = list(x = as.Date("2020-03-28"), y = 3.8, xref = "x", yref = "y", 
                                    text = "Limit of Detection", showarrow = FALSE))

        

        p1
## Warning: `arrange_()` is deprecated as of dplyr 0.7.0.
## Please use `arrange()` instead.
## See vignette('programming') for more help
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.
## Warning: Ignoring 1 observations
        p2
## Warning: `group_by_()` is deprecated as of dplyr 0.7.0.
## Please use `group_by()` instead.
## See vignette('programming') for more help
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.

Combine the two main plot pieces as a subplot

#seperate n1 and n2 frames by site
#n1
wrf_a_only_n1 <- subset(only_n1, Facility == "WRF A")
wrf_b_only_n1 <- subset(only_n1, Facility == "WRF B")
wrf_c_only_n1 <- subset(only_n1, Facility == "WRF C")

#n2
wrf_a_only_n2 <- subset(only_n2, Facility == "WRF A")
wrf_b_only_n2 <- subset(only_n2, Facility == "WRF B")
wrf_c_only_n2 <- subset(only_n2, Facility == "WRF C")


#rejoin the old data frames then seperate in to averages for each plant. 
wrfa_both <- full_join(wrf_a_only_n1, wrf_a_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))
## Joining, by = c("date", "cases_cum_clarke", "new_cases_clarke", "X7_day_ave_clarke", "cases_per_100000_clarke", "Facility", "collection_num", "target", "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "mean_total_copies", "sd_total_copies", "log_copy_per_L")
wrfb_both <- full_join(wrf_b_only_n1, wrf_b_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))
## Joining, by = c("date", "cases_cum_clarke", "new_cases_clarke", "X7_day_ave_clarke", "cases_per_100000_clarke", "Facility", "collection_num", "target", "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "mean_total_copies", "sd_total_copies", "log_copy_per_L")
wrfc_both <- full_join(wrf_c_only_n1, wrf_c_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))
## Joining, by = c("date", "cases_cum_clarke", "new_cases_clarke", "X7_day_ave_clarke", "cases_per_100000_clarke", "Facility", "collection_num", "target", "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "mean_total_copies", "sd_total_copies", "log_copy_per_L")
#get max date
maxdate <- max(wrfa_both$date)
mindate <- min(wrfa_both$date)

Build loess smoothing figures figures

This makes the individual plots

#**************************************WRF A PLOT**********************************************
#add trendlines 
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_botha <- ggplot(wrfa_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_botha<<-..y..), method = "loess", color = '#1B9E77', 
              span = 0.6, n = 198)
## Warning: Ignoring unknown aesthetics: outfit
#look at the fits to align dates and total observations
#both
extract_botha
## `geom_smooth()` using formula 'y ~ x'

fit_botha
##   [1] 12.73880 12.75184 12.76454 12.77692 12.78901 12.80085 12.81246 12.82387
##   [9] 12.83511 12.84621 12.85719 12.86809 12.87892 12.88973 12.90054 12.91129
##  [17] 12.92191 12.93238 12.94268 12.95278 12.96268 12.97235 12.98178 12.99094
##  [25] 12.99983 13.00841 13.01668 13.02461 13.03219 13.03954 13.04678 13.05390
##  [33] 13.06088 13.06771 13.07438 13.08087 13.08716 13.09324 13.09909 13.10471
##  [41] 13.11007 13.11517 13.11998 13.12449 13.12869 13.13256 13.13609 13.13926
##  [49] 13.14205 13.14447 13.14647 13.14807 13.14923 13.14994 13.15019 13.14997
##  [57] 13.14926 13.14783 13.14550 13.14233 13.13840 13.13375 13.12846 13.12258
##  [65] 13.11618 13.10932 13.10207 13.09447 13.08661 13.07853 13.07030 13.06199
##  [73] 13.05366 13.04536 13.03717 13.02914 13.02133 13.01381 13.00462 12.99209
##  [81] 12.97673 12.95908 12.93964 12.91893 12.89746 12.87577 12.85435 12.83374
##  [89] 12.81445 12.79699 12.78188 12.76964 12.75958 12.75065 12.74280 12.73601
##  [97] 12.73022 12.72540 12.72152 12.71853 12.71640 12.71509 12.71457 12.71478
## [105] 12.71571 12.71730 12.71953 12.72235 12.72572 12.72961 12.73398 12.73880
## [113] 12.74401 12.74959 12.75550 12.76170 12.76816 12.77482 12.78167 12.78865
## [121] 12.79360 12.79473 12.79258 12.78768 12.78056 12.77176 12.76181 12.75125
## [129] 12.74060 12.73040 12.72118 12.71349 12.70784 12.70478 12.70483 12.70401
## [137] 12.69914 12.69198 12.68432 12.67790 12.67452 12.67592 12.68187 12.69071
## [145] 12.70208 12.71562 12.73097 12.74777 12.76568 12.78432 12.80335 12.82241
## [153] 12.84113 12.85916 12.87614 12.89172 12.90693 12.92297 12.93970 12.95702
## [161] 12.97478 12.99288 13.01118 13.02956 13.04790 13.06607 13.08394 13.10140
## [169] 13.11904 13.13738 13.15621 13.17532 13.19448 13.21350 13.23214 13.25060
## [177] 13.26920 13.28794 13.30680 13.32580 13.34491 13.36414 13.38348 13.40293
## [185] 13.42248 13.44213 13.46187 13.48169 13.50160 13.52159 13.54165 13.56177
## [193] 13.58196 13.60221 13.62251 13.64285 13.66325 13.68368
#assign fits to a vector
both_trenda <- fit_botha

#extract y min and max for each
limits_botha <- ggplot_build(extract_botha)$data
## `geom_smooth()` using formula 'y ~ x'
limits_botha <- as.data.frame(limits_botha)
both_ymina <- limits_botha$ymin
both_ymaxa <- limits_botha$ymax

#reassign dataframes (just to be safe)
work_botha <- wrfa_both

#fill in missing dates to smooth fits
work_botha <- work_botha %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_botha <- work_botha$date

#create a new smooth dataframe to layer
smooth_frame_botha <- data.frame(date_vec_botha, both_trenda, both_ymina, both_ymaxa)
#WRF A
#plot smooth frames
p_wrf_a <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_botha, y = ~both_trenda,
                    data = smooth_frame_botha,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_botha,
                                  '</br> Median Log Copies: ', round(both_trenda, digits = 2)),
                    line = list(color = '#1B9E77', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_botha, ymin = ~both_ymina, ymax = ~both_ymaxa,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_botha, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxa, digits = 2),
                                  '</br> Min Log Copies: ', round(both_ymina, digits = 2)),
                    name = "",
                    fillcolor = '#1B9E77',
                    line = list(color = '#1B9E77')) %>%
                layout(yaxis = list(title = "Total Log SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF A") %>%
    plotly::add_segments(x = as.Date("2020-06-24"), 
                                          xend = as.Date("2020-06-24"), 
                                          y = ~min(both_ymina), yend = ~max(both_ymaxa),
                                          opacity = 0.35,
                                          name = "Bars Repoen",
                                          hoverinfo = "text",
                                          text = "</br> Bars Reopen",
                                                 "</br> 2020-06-24",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
    plotly::add_segments(x = as.Date("2020-07-09"), 
                                          xend = as.Date("2020-07-09"), 
                                          y = ~min(both_ymina), yend = ~max(both_ymaxa),
                                          opacity = 0.35,
                                          name = "Mask Mandate",
                                          hoverinfo = "text",
                                          text = "</br> Mask Mandate",
                                                 "</br> 2020-07-09",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
    plotly::add_segments(x = as.Date("2020-08-20"), 
                                          xend = as.Date("2020-08-20"), 
                                          y = ~min(both_ymina), yend = ~max(both_ymaxa),
                                          opacity = 0.35,
                                          name = "</br> Classes Begin",
                                                 "</br> 2020-08-20",
                                          hoverinfo = "text",
                                          text = "Classes Begin",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
        plotly::add_segments(x = as.Date("2020-10-03"), 
                                          xend = as.Date("2020-10-03"), 
                                          y = ~min(both_ymina), yend = ~max(both_ymaxa),
                                          opacity = 0.35,
                                          name = "</br> First Home Football Game",
                                                 "</br> 2020-10-03",
                                          hoverinfo = "text",
                                          text = "First Home Football Game",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfa_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#1B9E77', size = 6, opacity = 0.65))

p_wrf_a
save(p_wrf_a, file = "./plotly_objs/p_wrf_a.rda")
#**************************************WRF B PLOT**********************************************
#add trendlines 
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_bothb <- ggplot(wrfb_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_bothb<<-..y..), method = "loess", color = '#D95F02', 
              span = 0.6, n = 198)
## Warning: Ignoring unknown aesthetics: outfit
#look at the fits to align dates and total observations
#both
extract_bothb
## `geom_smooth()` using formula 'y ~ x'

fit_bothb
##   [1] 12.69081 12.67514 12.65961 12.64430 12.62927 12.61462 12.60041 12.58674
##   [9] 12.57367 12.56128 12.54966 12.53888 12.52902 12.52016 12.51237 12.50559
##  [17] 12.49963 12.49444 12.48996 12.48611 12.48285 12.48010 12.47781 12.47592
##  [25] 12.47437 12.47308 12.47201 12.47109 12.47026 12.46997 12.47069 12.47234
##  [33] 12.47484 12.47811 12.48209 12.48669 12.49185 12.49748 12.50351 12.50986
##  [41] 12.51647 12.52325 12.53013 12.53703 12.54389 12.55061 12.55714 12.56339
##  [49] 12.56929 12.57476 12.58160 12.59136 12.60359 12.61785 12.63370 12.65068
##  [57] 12.66835 12.68627 12.70399 12.72106 12.73705 12.75150 12.76397 12.77402
##  [65] 12.78418 12.79705 12.81232 12.82965 12.84871 12.86920 12.89077 12.91310
##  [73] 12.93586 12.95873 12.98139 13.00350 13.02474 13.04478 13.06331 13.07998
##  [81] 13.09448 13.10648 13.11565 13.12166 13.12420 13.12077 13.10982 13.09244
##  [89] 13.06969 13.04266 13.01243 12.98006 12.94665 12.91327 12.88099 12.85090
##  [97] 12.82406 12.80157 12.78450 12.76947 12.75253 12.73390 12.71378 12.69241
## [105] 12.66998 12.64671 12.62283 12.59853 12.57404 12.54957 12.52534 12.50155
## [113] 12.47843 12.45619 12.43504 12.41519 12.39687 12.38028 12.36564 12.35316
## [121] 12.34306 12.33556 12.33085 12.32917 12.33073 12.33573 12.34618 12.36294
## [129] 12.38452 12.40947 12.43629 12.46352 12.48968 12.51328 12.53286 12.55063
## [137] 12.56957 12.58935 12.60964 12.63009 12.65038 12.67018 12.69046 12.71226
## [145] 12.73538 12.75957 12.78463 12.81033 12.83645 12.86277 12.88905 12.91509
## [153] 12.94065 12.96552 12.98948 13.01229 13.03497 13.05849 13.08263 13.10714
## [161] 13.13180 13.15635 13.18057 13.20422 13.22706 13.24886 13.26937 13.28837
## [169] 13.30683 13.32567 13.34459 13.36327 13.38141 13.39869 13.41479 13.43000
## [177] 13.44478 13.45913 13.47305 13.48652 13.49956 13.51214 13.52427 13.53594
## [185] 13.54715 13.55790 13.56817 13.57797 13.58730 13.59613 13.60448 13.61234
## [193] 13.61970 13.62655 13.63291 13.63875 13.64407 13.64888
#assign fits to a vector
both_trendb <- fit_bothb

#extract y min and max for each
limits_bothb <- ggplot_build(extract_bothb)$data
## `geom_smooth()` using formula 'y ~ x'
limits_bothb <- as.data.frame(limits_bothb)
both_yminb <- limits_bothb$ymin
both_ymaxb <- limits_bothb$ymax

#reassign dataframes (just to be safe)
work_bothb <- wrfb_both

#fill in missing dates to smooth fits
work_bothb <- work_bothb %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothb <- work_bothb$date

#create a new smooth dataframe to layer
smooth_frame_bothb <- data.frame(date_vec_bothb, both_trendb, both_yminb, both_ymaxb)
#WRF B
#plot smooth frames
p_wrf_b <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_bothb, y = ~both_trendb,
                    data = smooth_frame_bothb,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothb,
                                  '</br> Median Log Copies: ', round(both_trendb, digits = 2)),
                    line = list(color = '#D95F02', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothb, ymin = ~both_yminb, ymax = ~both_ymaxb,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothb, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxb, digits = 2),
                                  '</br> Min Log Copies: ', round(both_yminb, digits = 2)),
                    name = "",
                    fillcolor = '#D95F02',
                    line = list(color = '#D95F02')) %>%
                layout(yaxis = list(title = "Total Log SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF B") %>%
    plotly::add_segments(x = as.Date("2020-06-24"), 
                                          xend = as.Date("2020-06-24"), 
                                          y = ~min(both_yminb), yend = ~max(both_ymaxb),
                                          opacity = 0.35,
                                          name = "Bars Repoen",
                                          hoverinfo = "text",
                                          text = "</br> Bars Reopen",
                                                 "</br> 2020-06-24",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
    plotly::add_segments(x = as.Date("2020-07-09"), 
                                          xend = as.Date("2020-07-09"), 
                                          y = ~min(both_yminb), yend = ~max(both_ymaxb),
                                          opacity = 0.35,
                                          name = "Mask Mandate",
                                          hoverinfo = "text",
                                          text = "</br> Mask Mandate",
                                                 "</br> 2020-07-09",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
    plotly::add_segments(x = as.Date("2020-08-20"), 
                                          xend = as.Date("2020-08-20"), 
                                          y = ~min(both_yminb), yend = ~max(both_ymaxb),
                                          opacity = 0.35,
                                          name = "</br> Classes Begin",
                                                 "</br> 2020-08-20",
                                          hoverinfo = "text",
                                          text = "Classes Begin",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
        plotly::add_segments(x = as.Date("2020-10-03"), 
                                          xend = as.Date("2020-10-03"), 
                                          y = ~min(both_yminb), yend = ~max(both_ymaxb),
                                          opacity = 0.35,
                                          name = "</br> First Home Football Game",
                                                 "</br> 2020-10-03",
                                          hoverinfo = "text",
                                          text = "First Home Football Game",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfb_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#D95F02', size = 6, opacity = 0.65))

p_wrf_b
save(p_wrf_b, file = "./plotly_objs/p_wrf_b.rda")

#**************************************WRF C PLOT********************************************** #add trendlines #extract data from geom_smooth # *********************************span 0.6*********************************** #*****************Must always update the n = TOTAL NUMBER OF DAYS*************************

extract_bothc <- ggplot(wrfc_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_bothc<<-..y..), method = "loess", color = '#E7298A', 
              span = 0.6, n = 198)
## Warning: Ignoring unknown aesthetics: outfit
#look at the fits to align dates and total observations
#both
extract_bothc
## `geom_smooth()` using formula 'y ~ x'

fit_bothc
##   [1] 11.15689 11.20554 11.25326 11.30005 11.34591 11.39083 11.43482 11.47787
##   [9] 11.51999 11.56117 11.60141 11.64071 11.67908 11.71649 11.75297 11.78850
##  [17] 11.82308 11.85672 11.88942 11.92116 11.95195 11.98179 12.01073 12.03880
##  [25] 12.06598 12.09227 12.11763 12.14206 12.16553 12.18804 12.20956 12.23007
##  [33] 12.24956 12.26802 12.28542 12.30175 12.31686 12.33064 12.34312 12.35434
##  [41] 12.36432 12.37311 12.38073 12.38722 12.39262 12.39694 12.40024 12.40253
##  [49] 12.40387 12.40426 12.40376 12.40239 12.40019 12.39719 12.39343 12.38893
##  [57] 12.38372 12.37786 12.37135 12.36425 12.35658 12.34838 12.33967 12.33050
##  [65] 12.31931 12.30475 12.28712 12.26673 12.24389 12.21891 12.19209 12.16374
##  [73] 12.13417 12.10368 12.07259 12.04121 12.00983 11.97876 11.94832 11.91882
##  [81] 11.89055 11.86383 11.83897 11.81627 11.79603 11.77348 11.74439 11.70989
##  [89] 11.67113 11.62925 11.58540 11.54071 11.49634 11.45343 11.41312 11.37654
##  [97] 11.34486 11.31920 11.30072 11.28745 11.27662 11.26814 11.26194 11.25793
## [105] 11.25604 11.25617 11.25825 11.26221 11.26794 11.27538 11.28445 11.29506
## [113] 11.30713 11.32057 11.33532 11.35128 11.36837 11.38652 11.40564 11.42565
## [121] 11.44647 11.46802 11.49021 11.51297 11.53621 11.55986 11.58459 11.61092
## [129] 11.63851 11.66701 11.69607 11.72536 11.75452 11.78322 11.81110 11.83737
## [137] 11.86195 11.88562 11.90912 11.93323 11.95870 11.98631 12.01510 12.04370
## [145] 12.07218 12.10064 12.12915 12.15781 12.18669 12.21588 12.24547 12.27554
## [153] 12.30617 12.33745 12.36947 12.40231 12.43687 12.47360 12.51183 12.55092
## [161] 12.59023 12.62910 12.66690 12.70296 12.73666 12.76733 12.79433 12.81702
## [169] 12.83744 12.85781 12.87773 12.89682 12.91468 12.93093 12.94516 12.95774
## [177] 12.96929 12.97978 12.98920 12.99753 13.00476 13.01087 13.01585 13.01968
## [185] 13.02234 13.02381 13.02409 13.02316 13.02099 13.01757 13.01290 13.00694
## [193] 12.99969 12.99113 12.98124 12.97001 12.95742 12.94346
#assign fits to a vector
both_trendc <- fit_bothc

#extract y min and max for each
limits_bothc <- ggplot_build(extract_bothc)$data
## `geom_smooth()` using formula 'y ~ x'
limits_bothc <- as.data.frame(limits_bothc)
both_yminc <- limits_bothc$ymin
both_ymaxc <- limits_bothc$ymax

#reassign dataframes (just to be safe)
work_bothc <- wrfc_both

#fill in missing dates to smooth fits
work_bothc <- work_bothc %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothc <- work_bothc$date

#create a new smooth dataframe to layer
smooth_frame_bothc <- data.frame(date_vec_bothc, both_trendc, both_yminc, both_ymaxc)
#WRF C
#plot smooth frames
p_wrf_c <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_bothc, y = ~both_trendc,
                    data = smooth_frame_bothc,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothc,
                                  '</br> Median Log Copies: ', round(both_trendc, digits = 2)),
                    line = list(color = '#E7298A', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothc, ymin = ~both_yminc, ymax = ~both_ymaxc,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothc, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxc, digits = 2),
                                  '</br> Min Log Copies: ', round(both_yminc, digits = 2)),
                    name = "",
                    fillcolor = '#E7298A',
                    line = list(color = '#E7298A')) %>%
                layout(yaxis = list(title = "Total Log SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF C") %>%
    plotly::add_segments(x = as.Date("2020-06-24"), 
                                          xend = as.Date("2020-06-24"), 
                                          y = ~min(both_yminc), yend = ~max(both_ymaxc),
                                          opacity = 0.35,
                                          name = "Bars Repoen",
                                          hoverinfo = "text",
                                          text = "</br> Bars Reopen",
                                                 "</br> 2020-06-24",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
    plotly::add_segments(x = as.Date("2020-07-09"), 
                                          xend = as.Date("2020-07-09"), 
                                          y = ~min(both_yminc), yend = ~max(both_ymaxc),
                                          opacity = 0.35,
                                          name = "Mask Mandate",
                                          hoverinfo = "text",
                                          text = "</br> Mask Mandate",
                                                 "</br> 2020-07-09",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
    plotly::add_segments(x = as.Date("2020-08-20"), 
                                          xend = as.Date("2020-08-20"), 
                                          y = ~min(both_yminc), yend = ~max(both_ymaxc),
                                          opacity = 0.35,
                                          name = "</br> Classes Begin",
                                                 "</br> 2020-08-20",
                                          hoverinfo = "text",
                                          text = "Classes Begin",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
        plotly::add_segments(x = as.Date("2020-10-03"), 
                                          xend = as.Date("2020-10-03"), 
                                          y = ~min(both_yminc), yend = ~max(both_ymaxc),
                                          opacity = 0.35,
                                          name = "</br> First Home Football Game",
                                                 "</br> 2020-10-03",
                                          hoverinfo = "text",
                                          text = "First Home Football Game",
                                          showlegend = FALSE,
                                          line = list(color = "black", dash = "dash")) %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfc_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#E7298A', size = 6, opacity = 0.65))

p_wrf_c
save(p_wrf_c, file = "./plotly_objs/p_wrf_c.rda")
save(wrfa_both, file = "./plotly_objs/wrfa_both.rda")
save(wrfb_both, file = "./plotly_objs/wrfb_both.rda")
save(wrfc_both, file = "./plotly_objs/wrfc_both.rda")
save(date_vec_botha, file = "./plotly_objs/date_vec_botha.rda")
save(date_vec_bothb, file = "./plotly_objs/date_vec_bothb.rda")
save(date_vec_bothc, file = "./plotly_objs/date_vec_bothc.rda")
save(both_ymina, file = "./plotly_objs/both_ymina.rda")
save(both_ymaxa, file = "./plotly_objs/both_ymaxa.rda")

save(both_yminb, file = "./plotly_objs/both_yminb.rda")
save(both_ymaxb, file = "./plotly_objs/both_ymaxb.rda")

save(both_yminc, file = "./plotly_objs/both_yminc.rda")
save(both_ymaxc, file = "./plotly_objs/both_ymaxc.rda")